Living polymers in a size-asymmetric electrolyte 
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A living polymers transition is found in molecular dynamics simulations of a charge-symmetric 
size-asymmetric electrolyte with no anisotropic interactions. The fluid has strong polymeric char- 
acter at low temperatures, where it consists of large, alternating-charge linear chains and rings 
in chemical equilibrium. A mean-field theory of chain association is consistent with the existence 
of such a transition. In the polymeric phase the system is very weakly conducting or electrically 
insulating. 

PACS numbers: 61.20.Qg, 66.10.Ed, 61.25.Hq 



Ionic systems are rather ubiquitous in nature, ranging 
from molten salts to the electrolyte solutions crucial for 
the functioning of living cells, to even pure water at very 
high densities jjj. A decisive step in understanding how 
the long-range character of the Coulomb potential deter- 
mines the behavior of ionic systems has long proved to 
be the Debye-Hiickel theory of screening by counterions 
H . The introduction of the physically appealing notion 
of ionic pairing led to an additional refinement of these 
ideas and has allowed the recent development of success- 
ful theories for criticality in simple electrolyte models || . 
The formation of ionic pairs and higher order clusters was 
initially observed and characterized in symmetric elec- 
trolytes j| and it is now understood to be an essential 
feature of ionic systems || |[ 0, |[ [||. Recent work has 
revealed that ionic clustering is even more pronounced in 
size-asymmetric model electrolytes and it is associated 
with an unexpected behavior of their phase coexistence 
regions ffl]. Size-asymmetry is of course the norm in 
numerous circumstances, e.g. colloid science pf]| or sys- 
tems with significant ionic character such as LiCl, NaH 
etc. We present molecular dynamics simulations of a size- 
asymmetric electrolyte with isotropic interactions, where 
strong ionic clustering produces a continuous polymeriza- 
tion transition to a phase with very weakly conducting 
or insulating character. 

We study a system of N/2 positive charges +q and 
N/2 negative charges — q interacting through potentials 



Vij(r) = Uij(r) + 



(1) 



(i,j = +,— ). Dq is the vacuum dielectric permeability, 
and Uij{r) are exponential- 6 potentials, 



u ij( r ) = e 



Ae 



D 



(2) 



A 



with well-depth e corresponding to distance 
6e a /(a - 6), B = a/(a - 6) and a = 13. The Born- 
Mayer exponential term provides a good description of 
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the repulsive short range interactions that occur between 
closed-shell atoms or ions, e.g. in molten salts When 
modified to include many-body effects such simple ionic 
models appear to have a rather wide range of applicabil- 
ity |Q . We set q = 1 in units of the electron charge e and 
impose an extreme size-asymmetry by choosing 7' ++ = 

and = r /2, i.e. the positive ions are point 

charges. We choose as relevant length scale the position 

of the minimum of (r), which we denote by , 

d_l ~ 0.67r^ This is appropriate at low temperatures 

where it should allow a comparison with the charged hard 
core models usually studied. The total reduced number 
density is p* — pd+_, p = N/V, and the reduced tem- 
perature T* = k B T/E , with E = q 2 /4TrD d + -. Most 
of the simulations were performed with N — 2048 par- 
ticles in the microcanonical (NVE) ensemble, in a box 
with periodic boundary conditions, at reduced density 
p* = 0.12 and reduced temperatures T* in the range 
0.0065 to 0.0365. The reduced simulation box length was 
L* = L/d+- ~ 26; simulations with L* ~ 35 (N = 5324) 
were also performed to verify the influence of finite size 
effects. The energy scale of the short range potentials 
was set to c/Eq ~ 2.4 x 10~ 3 . The Coulomb interac- 
tions were handled using the Ewald summation technique 
with conducting boundary conditions. The system was 
typically equilibrated for 5 x 10 4 — 10 5 time steps, using 
successively velocity rescaling, Andersen thermostat and 
normal (NVE) molecular dynamics. The accumulation 
runs consisted of at least 5 x 10 5 time steps. 

At the lowest temperatures the pair correlation func- 
tions signal a strongly structured fluid, Fig. 1. The 

opposite-charge pair correlation function (r) exhibits 

a very sharp peak at small distances, one that might 
be expected for example if significant ionic pairing takes 
place. An analysis of the coordination number as a func- 
tion of distance, N + -(r) = [ft/2) g + -(r/)4:Trr/ 2 dr!, Fig. 
1 (inset), reveals however the presence of two opposite- 
charge nearest neighbors, an indication of alternating- 
charge chain association. The importance of such 
polymeric structures for the thermodynamics of elec- 
trolytes has been suggested before in the context of size- 
symmetric systems [gj and has been recently reinforced 
by Monte Carlo simulations of size-asymmetric hard- 
core electrolytes |6[ 0. Due to the large size-asymmetry 
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FIG. 1: Opposite-charge pair correlation functions g+-(r) at 
T* = 0.0065 (dotted line), T* = 0.016 (solid line), T* = 
0.024 (dashed line) and corresponding coordination numbers 
N+-(r) (inset). 



their contribution is dominant here, apparently leading 
- see Fig. 1 - to a continuous transition between a non- 
structured isotropic fluid at high temperatures to a struc- 
tured isotropic fluid at low temperatures. 

The formation of large structures in ionic systems can 
be understood using appropriate clustering definitions as 
first introduced for ionic systems by Gillan [Q. Gillan's 
original analysis of cluster populations labels an ion as 
belonging to a cluster if the smallest distance between 
the ion and other cluster ions is less than a certain dis- 
tance Re- We use Rc — 1.5<i + _, which is a good approx- 
imation for the position of the first minimum of (?*) . 

We note however that the low temperatures analysis is 
insensitive to the choice of Rc in a significant range. 
The same is true if we adopt a more restrictive defini- 
tion of clustering, i.e. one that requires that an ion be 
'linked' to a cluster through an ion of opposite charge. 
We find that at low temperatures the system consists of 
large alternating-charge chains and rings that behave as 
living polymers, i.e. break and recombine continuously 
leading to a temperature-dependent equilibrium size dis- 
tribution. The lifetime of these linear polymers, particu- 
larly the smaller ones, is rather long, which explains the 
extremely long times needed to equilibrate the system. 

Fig. 2 shows snapshots of the motion of such a small 
ring at the lowest temperature simulated, which diffuses 
as a well defined entity and exhibits the vibrational and 
rotational motions specific of a regular polymer. 

The association of neutral monomers in linear poly- 
mers under conditions of chemical equilibrium leads to 
a continuous polymerization transition accompanied by 
a maximum of the specific heat signaling the onset of 
chain formation Jl3[ . The position of the maximum is 
a measure of the association, or 'bond' energy within a 
chain. Such a transition is also apparently taking place in 




FIG. 2: Snapshots of the motion of a small alternating-charge 
ring at intervals of 150000 time steps; T* = 0.0065. 
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FIG. 3: Specific heat cv from simulations (circles) and mean- 
field theory (solid line) after a temperature rescaling (see 
text). Inset: average chain size n a from simulations with 
N = 2048 particles (circles), N = 5324 (filled diamonds) and 
mean-field theory (solid line); arrow indicates T£ . 



the ionic size- asymmetric system that we are simulating. 
We therefore calculate the constant volume specific heat, 
cy = (dE/dT)y/N, as a function of temperature in the 
microcanonical ensemble using the equilibrium kinetic 
energy fluctuations , (/C 2 ) - (/C) 2 = (3/2)iVfe|T 2 (l - 
3&b/2cv) Q. The maximum of cv, Fig. 3, places the 
continuous polymerization transition at T* ~ 0.016, in 
agreement with the position of energy E = E(T) and 
pressure p — p(T) inflection points. It is remarkable 
that in the living polymers systems previously studied 
the association of monomers into chains is driven by 
strong anisotropic interactions, which many times are of 
a dipolar nature fllf . While in the present case no such 
anisotropy is present, isotropic long range (Coulomb) in- 
teractions mediated by extreme size- asymmetry lead to 
similar characteristic behavior. This could be perhaps 
attributed to an enhanced tendency for dipolar pairing 
due to size-asymmetry. However, although the chain 
size distribution above the polymerization temperature 
is peaked at n — 2 (dipolar pairs), the distribution of ion 
fractions making up identical size chains is not, with the 
value for pairs remaining very small, of the order 10~ 2 . 
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Nevertheless, energy considerations also suggest that the 
formation of chains is favored in systems with large size- 
asymmetry ^,0). If the propensity for chain formation 
is assumed to prevail for such electrolytes, a mean-field 
theory for charged associating chains can be formulated, 
as in the case of neutral living polymers [fl6| . At the sim- 
plest level such a theory (see below) does not include the 
size- asymmetry per se, but only implicitly in the chain 
association assumption. The system is assumed to be 
a mixture of alternating-charge chains which are electri- 
cally neutral when containing equal numbers of positive 
and negative ions and with charge ±g otherwise. At low 
densities the free energy of such a system can be written 
as 

F 3 
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where p\ n and p 2n are densities of uncharged and charged 
chains, respectively. K n is the partition function of an 
n— chain and p c is the total density of charged chains 
Pc = J2„= (P2n+i + PiTn+i)- The Debye-Hiickel term 
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includes an effective hard core d, that we set to d^ in- 
dependent of polymerization effects; k 2 = q 2 p c j D^ksT '. 
For the n— chain partition function we use K n = 
[K2(T)] n ~ 1 , with K2{T) the dimer association constant 
in the low temperature Bjerrum approximation, KziT) — 
A-Kd\_T*exp(l/T*) ||, neglecting interactions beyond 
nearest neighbors along a chain. (Here we treat the 

short range potential (r) as effectively hard core, but 

this should introduce only small errors at low tempera- 
tures.) We note that if only H — 'dimerization' is allowed, 
the above theory reduces to the so-called Debye-Huckel- 
Bjerrum (DHBj) theory ||. The energy of the system is 
the sum of kinetic, 'bond' and ionic terms, 

E 3 

e= — = -k B T+(p- p t )e 2 + e DH (p c ) (5) 

with p t the total chains density, pt = p n , e 2 = 

—E + ksT the average H — 'bond' energy and eoH the 
Debye-Hiickel contribution, 

E T* M+_) 3 



eDfi(pc) = - 
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The chain size distribution p n is then obtained by mini- 
mizing / at fixed density p and temperature T, with the 
constraint of particle conservation, X^n=i n P n = P> an< ^ 
charge neutrality, £„ =0 [ftLfi _ P2 n +l] = °- We obtain 
p° 2n = K^e 2nX (7) 

PUi = p *n + x = K^e^exp (-^) (8) 



where pj n+1 = p^n+i ^ s a consequence of charge neutral- 
ity and A can be interpreted as the monomer chemical 
potential. The particle conservation condition can be re- 
cast as an equation for the charge density p Cl which after 
being solved numerically also allows the calculation of pt ■ 
Substitution into the energy equation yields the specific 
heat through pcy = de/dT\ p . 

We find that cy displays the characteristic maximum 
]l3[ of a living polymers transition at T* ~ 0.11. This is 
higher than the transition temperature obtained in the 
simulations, but it is not surprising given the many ap- 
proximations involved. In particular, the K n approxi- 
mation used leads to a large overestimate even for the 
energy of long, straight alternating-charge chains, while 
size-asymmetry dependent packing and solvation effects, 
neglected here, should decrease the polymerization tem- 
perature by weakening the H — 'bonds' [0. After a 
n+1 suitable rescaling of the temperature the calculated cy 
is qualitatively compared with the simulation results in 
Fig/.g-jB. We note that the height of the cy maximum is 
very well reproduced, along with the higher temperatures 
behavior. The low temperature deviations however are 
significant, probably due to the neglect of ring structures 
that appear to be very important below the transition 
temperature. In fact the cluster analysis yields mostly 
neutral structures at low temperatures, which is not well 
reproduced by the above mean-field model. 

Another signature of a continuous polymerization tran- 
sition is the rapid chain growth below T* [[13). We 
show in Fig. 3 (inset) the average chain size n a — 
S n =i n Pnl y~l n — l P n f° un d in the simulations and with 
the above mean-field treatment after the same rescaling 
of the temperature used for cy . The chain size at T* is in 
good agreement with the theory, as it is the overall behav- 
ior. The slower rise of n a in the simulations than in the 
theory at low temperatures may be due to a depletion of 
linear chains due to rings formation, as noted before. It is 
worth pointing out that the finite size effects that appear 
to affect n a for the smaller system size occur only at the 
lowest temperatures, far away from T* , We also estimate 
an effective lifetime r of the polymeric chains by assum- 
ing an initial exponential decay exp(—t/r) of the time 
correlation function C(At) — (pt(to)pt(to + At)) — (pt) 2 
of the total number of chains p t ; r is for example of the 
order of 20 MD time steps at T* = 0.01 (compared with 
~ 5 x 10 5 total simulation steps), although it may be 
significantly longer for the smaller chains and rings. 

The formation of large chains and rings and the dom- 
inance of the neutral ones has important consequences 
for the conduction properties of the polymeric phase. We 
calculate the electrical conductivity a using the Green- 
Kubo expression 



a(t) = 



1 



3Vk B T 



(j(*')-j(0)>df 



(9) 



where j(t) = X)fcli9fc v fcW i s the charge current; a 
is evaluated from the large times value of a(t), a — 
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FIG. 4: Electrical conductivity for T* > T* (circles) and fit 
with a — aoexp(—E^/T*) (dotted line). 



limt_»oo Below the transition temperature a(t) ex- 
hibits very strong oscillatory behavior and decays ap- 
parently to zero, i.e. the polymeric phase is very weakly 
conducting or insulating; the conductivities above T* are 
shown in Fig. 4. As the temperature increases the con- 
ductivity also increases exhibiting an activated behavior 
well described by a = a a exp{~E* b /T*), with E£ ~ 0.1. 
The conductivity increase is driven by an increasing 
charged clusters density which occurs through the break- 
up of H — 'bonds'. Ef, can therefore be interpreted as an 
effective 'bond' strength, yielding T*/E£ ~ 0.16, which 
is comparable with the value found in neutral living poly- 
mers of 0.25 |l3| . The possibility of an electrolyte insulat- 
ing phase has been discussed before in the light of strong 
ionic pairing at low densities, in the vapor region of sym- 



metric ionic systems p], [l9| . An insulating phase 
may arise here in connection with a living polymers tran- 
sition in the liquid region of the phase diagram. It should 
be noted that the very slow dynamics of the system at 
low temperatures, leading to the observed Arrhenius tem- 
perature dependence of the conductivity, may render the 
question of the conducting or insulating character of the 
polymeric phase rather moot. 

The chain associating tendency of size-asymmetric 
electrolytes, underscored by the continuous polymeriza- 
tion transition reported here, may interfere with the 
vapor-liquid coexistence and perhaps explain the ob- 
served reduction of the vapor-liquid phase separation re- 
gion in electrolytes as the size-asymmetry increases |j, [?| . 
The competition between chain association and liquid 
condensation has been studied before for the case of neu- 
tral monomers, where it was found that vapor-liquid co- 
existence can become metastable under strong chain as- 
sociation conditions p0[ . This is also worth exploring 
in the case of ionic systems, particularly because simu- 
lation results suggest that vapor-liquid coexistence may 
persist in charge-symmetric size-asymmetric electrolytes 
even for the extreme case of point cations || , while theo- 
retical calculations find no such coexistence jl0| . Finally, 
we note that polymeric behavior has been observed in 
molten BeC^ pi[ , where it was found to be due to the 
strong polarizability of the anion, coupled with the size- 
asymmetry of the ions. Wc find here that size- asymmetry 
alone may lead to living polymers behavior. 
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